home *** CD-ROM | disk | FTP | other *** search
/ Acorn RISC PD-CD 1 / Acorn RISC PD-CD 1.iso / languages / perl / lib / bigrat < prev    next >
Text File  |  1990-11-13  |  4KB  |  147 lines
  1. package bigrat;
  2. require "bigint";
  3.  
  4. # Arbitrary size rational math package
  5. #
  6. # Input values to these routines consist of strings of the form 
  7. #   m|^\s*[+-]?[\d\s]+(/[\d\s]+)?$|.
  8. # Examples:
  9. #   "+0/1"                          canonical zero value
  10. #   "3"                             canonical value "+3/1"
  11. #   "   -123/123 123"               canonical value "-1/1001"
  12. #   "123 456/7890"                  canonical value "+20576/1315"
  13. # Output values always include a sign and no leading zeros or
  14. #   white space.
  15. # This package makes use of the bigint package.
  16. # The string 'NaN' is used to represent the result when input arguments 
  17. #   that are not numbers, as well as the result of dividing by zero and
  18. #       the sqrt of a negative number.
  19. # Extreamly naive algorthims are used.
  20. #
  21. # Routines provided are:
  22. #
  23. #   rneg(RAT) return RAT                negation
  24. #   rabs(RAT) return RAT                absolute value
  25. #   rcmp(RAT,RAT) return CODE           compare numbers (undef,<0,=0,>0)
  26. #   radd(RAT,RAT) return RAT            addition
  27. #   rsub(RAT,RAT) return RAT            subtraction
  28. #   rmul(RAT,RAT) return RAT            multiplication
  29. #   rdiv(RAT,RAT) return RAT            division
  30. #   rmod(RAT) return (RAT,RAT)          integer and fractional parts
  31. #   rnorm(RAT) return RAT               normalization
  32. #   rsqrt(RAT, cycles) return RAT       square root
  33.  
  34. # Convert a number to the canonical string form m|^[+-]\d+/\d+|.
  35. sub main'rnorm { #(string) return rat_num
  36.     local($_) = @_;
  37.     s/\s+//g;
  38.     if (m#^([+-]?\d+)(/(\d*[1-9]0*))?$#) {
  39.     &norm($1, $3 ? $3 : '+1');
  40.     } else {
  41.     'NaN';
  42.     }
  43. }
  44.  
  45. # Normalize by reducing to lowest terms
  46. sub norm { #(bint, bint) return rat_num
  47.     local($num,$dom) = @_;
  48.     if ($num eq 'NaN') {
  49.     'NaN';
  50.     } elsif ($dom eq 'NaN') {
  51.     'NaN';
  52.     } elsif ($dom =~ /^[+-]?0+$/) {
  53.     'NaN';
  54.     } else {
  55.     local($gcd) = &'bgcd($num,$dom);
  56.     if ($gcd ne '+1') { 
  57.         $num = &'bdiv($num,$gcd);
  58.         $dom = &'bdiv($dom,$gcd);
  59.     } else {
  60.         $num = &'bnorm($num);
  61.         $dom = &'bnorm($dom);
  62.     }
  63.     substr($dom,0,1) = '';
  64.     "$num/$dom";
  65.     }
  66. }
  67.  
  68. # negation
  69. sub main'rneg { #(rat_num) return rat_num
  70.     local($_) = &'rnorm($_[0]);
  71.     tr/-+/+-/ if ($_ ne '+0/1');
  72.     $_;
  73. }
  74.  
  75. # absolute value
  76. sub main'rabs { #(rat_num) return $rat_num
  77.     local($_) = &'rnorm($_[0]);
  78.     substr($_,0,1) = '+';
  79.     $_;
  80. }
  81.  
  82. # multipication
  83. sub main'rmul { #(rat_num, rat_num) return rat_num
  84.     local($xn,$xd) = split('/',&'rnorm($_[0]));
  85.     local($yn,$yd) = split('/',&'rnorm($_[1]));
  86.     &norm(&'bmul($xn,$yn),&'bmul($xd,$yd));
  87. }
  88.  
  89. # division
  90. sub main'rdiv { #(rat_num, rat_num) return rat_num
  91.     local($xn,$xd) = split('/',&'rnorm($_[0]));
  92.     local($yn,$yd) = split('/',&'rnorm($_[1]));
  93.     &norm(&'bmul($xn,$yd),&'bmul($xd,$yn));
  94. }
  95.  
  96. # addition
  97. sub main'radd { #(rat_num, rat_num) return rat_num
  98.     local($xn,$xd) = split('/',&'rnorm($_[0]));
  99.     local($yn,$yd) = split('/',&'rnorm($_[1]));
  100.     &norm(&'badd(&'bmul($xn,$yd),&'bmul($yn,$xd)),&'bmul($xd,$yd));
  101. }
  102.  
  103. # subtraction
  104. sub main'rsub { #(rat_num, rat_num) return rat_num
  105.     local($xn,$xd) = split('/',&'rnorm($_[0]));
  106.     local($yn,$yd) = split('/',&'rnorm($_[1]));
  107.     &norm(&'bsub(&'bmul($xn,$yd),&'bmul($yn,$xd)),&'bmul($xd,$yd));
  108. }
  109.  
  110. # comparison
  111. sub main'rcmp { #(rat_num, rat_num) return cond_code
  112.     local($xn,$xd) = split('/',&'rnorm($_[0]));
  113.     local($yn,$yd) = split('/',&'rnorm($_[1]));
  114.     &bigint'cmp(&'bmul($xn,$yd),&'bmul($yn,$xd));
  115. }
  116.  
  117. # int and frac parts
  118. sub main'rmod { #(rat_num) return (rat_num,rat_num)
  119.     local($xn,$xd) = split('/',&'rnorm($_[0]));
  120.     local($i,$f) = &'bdiv($xn,$xd);
  121.     if (wantarray) {
  122.     ("$i/1", "$f/$xd");
  123.     } else {
  124.     "$i/1";
  125.     }   
  126. }
  127.  
  128. # square root by Newtons method.
  129. #   cycles specifies the number of iterations default: 5
  130. sub main'rsqrt { #(fnum_str[, cycles]) return fnum_str
  131.     local($x, $scale) = (&'rnorm($_[0]), $_[1]);
  132.     if ($x eq 'NaN') {
  133.     'NaN';
  134.     } elsif ($x =~ /^-/) {
  135.     'NaN';
  136.     } else {
  137.     local($gscale, $guess) = (0, '+1/1');
  138.     $scale = 5 if (!$scale);
  139.     while ($gscale++ < $scale) {
  140.         $guess = &'rmul(&'radd($guess,&'rdiv($x,$guess)),"+1/2");
  141.     }
  142.     "$guess";          # quotes necessary due to perl bug
  143.     }
  144. }
  145.  
  146. 1;
  147.